Model and Simulations of the Epitaxial Growth of Graphene on Non-Planar 6H-SiC 
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We study step flow growth of epitaxial graphene on 6H-SiC using a one dimensional kinetic 
Monte Carlo model. The model parameters are effective energy barriers for the nucleation and 
propagation of graphene at the SiC steps. When the model is applied to graphene growth on vicinal 
surfaces, a strip width distribution is used to characterize the surface morphology. Additional kinetic 
processes are included to study graphene growth on SiC nano-facets. Our main result is that the 
original nano-facet is fractured into several nano-facets during graphene growth. This phenomenon 
is characterized by the angle at which the fractured nano-facet is oriented with respect to the 
basal plane. The distribution of this angle across the surface is found to be related to the strip 
width distribution for vicinal surfaces. As the terrace propagation barrier decreases, the fracture 
angle distribution changes continously from two-sided Gaussian to one-sided power-law. Using this 
distribution, it will be possible to extract energy barriers from experiments and interpret the growth 
morphology quantitatively. 

PACS numbers: 68.35.-p, 68.55.-a, 81.15.Aa 



I. INTRODUCTION 

Epitaxial graphene grown by thermal decomposition of 
SiC is an attractive enabling technology for integrating 
graphene into silicon microelectronicsii^. However, scal- 
able production of high quality graphene film is still chal- 
lenging. Progress in this field requires understanding of 
the growth mechanism. Unfortunately, most experimen- 
tal results on epitaxial graphene growth are qualitative. 
The same is true for the interpretation of these results. 

First-principle approaches to this subject are difficult 
due to the complexity of the graphene/SiC interface^^^ 
and an unknown step structure. Previously, the present 
authors introduced a simple one dimensional kinetic 
Monte Carlo (KMC) model to study the epitaxial growth 
of graphene on vicinal surfaces of 6H-SiC^. Based on a 
study of graphene strip width distributions, the simula- 
tion results showed two distinct growth regimes domi- 
nated by "coalescence" processes and "climb-over" pro- 
cesses, respectively. 

Besides growth on SiC vicinal substrates, it has been 
shown recently that graphene grows spontaneously on 
nano-facetted SiC substrates^ . In this context, a 
"nano-facet" refers to a stable structure composed of 
several closely-spaced triple bilayer steps on the flat SiC 
surface. Suprisingly, perhaps, growth on a nano-facetted 
surface often leads to better quality graphene films on the 
adjacent terrace^' ^. Growth on nano-facets also plays an 
important role in graphene ribbon growth^. The purpose 
of this paper is to extend our previous study for vicinal 
substrates to include growth on nano-facetted substrates. 

Commercially available SiC substrates exhibt rough 
and scratched surfaces. An effective method to produce 
atomically flat vicinal surfaces is high temperature H2 
etching. For on-axis 6H-SiC substrates, this method pro- 
duces vicinal surfaces with triple bilayer steps. This type 
of microstep formation is related to etching kinetics and 



energy differences between different basal planes^^. This 
contrasts with off-axis 6H-SiC substrates where the H2 
etched surface often shows periodic structures of nano- 
facets. Depending on whether the substrate is tilted to- 
ward the (1100) or (1120) direction, the nano-facet makes 
an angle of ~ 25° or 13 — 14° from the basal plane^^^^. 
Besides H2 etching, controled nano-facets can be achieved 
by direct plasma etching of SiC surfaces^. 

Nano-facets can also form spontaneously when SiC is 
heated close to the graphitization temperature in non- 
vacuum environments. This type of step bunching is ob- 
served for both on- and off-axis SiC substrates. At higher 
growth pressure, a reduced silicon sublimation rate leads 
to a higher graphitization temperature. It is likely that 
at these elevated temperatures, the SiC vicinal surface 
of triple bilayer steps is further reconstructed into nano- 
facets to minimize the surface free energy. The forma- 
tion of graphene begins at these nano-facets where silicon 
atoms have fewer bonds^^^^. Given the previous discus- 
sion, SiC step bunching to produce nano-facets may oc- 
cur before graphene formation. This paper only considers 
the formation of graphene on pre-existing nano-facets re- 
gardless of their origin. We also neglect the difference 
between (llOn) and (112n) nano-facets. 

The structure of this paper is as follows. Section II is a 
brief review of our previous results for graphene growth 
on vicinal surfaces. In section III, we introduce our new 
model for graphene growth on a SiC surface with nano- 
facets. Section IV focuses on the distribution of fracture 
angles which can be compared directly with experiments 
to extract effective energy barriers. 



II. STRIP WIDTH DISTRIBUTION 

Our original model applies to graphene growth on a 
vicinal surface of 6II-SiC composed of triple bilayer steps. 
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FIG. 1. Graphene growth kinetic processes on a vicinal sur- 
face. 



barrier is small compared to the propagation barrier, i.e. 

< 0, all SiC steps decompose almost simultaneously 
at the beginning of growth, and the strip width distri- 
bution is Poisson (Fig. Efa)). In this case, coalescence 
processes occur in a short time interval right before the 
substrate is fully covered by graphene. Increasing the 
nucleation barrier relative to the propagation barrier re- 
duces the number of graphene strips on the surface and 
makes "climb-over" processes more dominant than co- 
alescences. This decreases the peak magnitude in the 
strip width distribution and shifts the peak position to 
the right in Fig. [Sl^b) and (c). Eventually, when AE is 
big enough, the step edges become indistinguishable as 
the graphene strip grows. The strip width distribution is 
then uniform(Fig. [2jd)). 



III. MODELING GROWTH ON A NANO-FACET 



From mass conservation, one graphene unit is obtained 
by decomposing one SiC step. Step flow growth from SiC 
triple bilayer steps produces strips of graphene parallel to 
the SiC step edges. This allows us to model the growth 
as one-dimensional. The detail of the atomic structure 
of the graphene/SiC interface does not play a role here, 
so we do not consider it in the model. Two model energy 
barriers are necessary: one is the barrier to graphene nu- 
cleation at a step edge (Fig. ma)-(b)). The other is the 
barrier to graphene strip propagation. Nucleation occurs 
at a rate rnuc = exp(— ^^nuc/^^)^ where ~ 10"*^^ 
s~-^ is an attempt frequency and T is the substrate tem- 
perature. Propagation occurs (Fig. [Hb)-(c)) at a rate 
^prop = I'o^^'Pi—Eprop/kT). The two barriers are "ef- 
fective" energy parameters which account for Si atoms 
sublimation, C atoms re-crystallization and subsequent 
graphene growth along the step edgei^^. 

After a series of propagation steps, two outcomes are 
possible for each graphene strip. One is that a growing 
graphene strip runs into another SiC step. We allow the 
growth to continue onto the adjacent upper terrace. This 
is called a "climb-over" process (Fig.[Hd)-(e)). The other 
outcome is that a growing graphene strip meets another 
graphene strip on the upper terrace. In this case, the 
two graphene strips connect to each other to form one 
continuous strip. We call this a "coalescence" process 
(Fig. [Hf)-(g)). Both climb-over and coalescence produce 
a step which is covered by a continuous graphene layer. A 
second graphene layer can nucleate at such a step under 
the first layer. We allow the subsequent growth of the 
second layer graphene to continue in the same way as 
the first layer. 

Our major finding from this model is the cross-over of 
the graphene strip width distribution shown in Fig. [2l 
This cross-over is the result of competition between the 
climb-over and coalescence processes. We define AE' = 
Enuc — ^prop and total coverage O = J^i^Q^^ where 0^ 
is the graphene coverage of layer i. When the nucleation 



In our modeling, a nano-facet is a group of triple bi- 
layer SiC steps with a fixed spacing between adjacent 
steps (Fig. [3j). Hence, the graphene growth proceeds 
also in step flow mode which can be modeled as one- 
dimensional. The nano-facet can have different angles 
with respect to the basal plane depending on both the 
orientation of the substrate and the graphene growth 
conditions. Adjusting the spacing width in the model al- 
lows different initial nano-facet angles. We assume that 
growth starts at the bottom of the nano-facet, converting 
one step into one graphene unit. The nano-facet nucle- 
ation process occurs at a rate rnuc = exp(— £^nuc/^2^) 
(Fig. [3l^a)-(b)). The nucleated graphene unit propa- 
gates upward along the nano-facet at a rate rprop = 
uoex.-p{—Epj-op/kT) (Fig. [3fb)-(c)). We keep the value of 
^prop = because incomplete graphene coverage on the 
nano-facet is rarely observed experimentally. In other 
words, growth on the nano-facet is very fast. 

The nano-facet nucleation and propagation processes 
defined here are very similar to those for vicinal surfaces. 
As soon as a nano-facet propagation event occurs, a sec- 
ond graphene layer can nucleate immediately under the 
first layer (Fig. [3ff)) and the growth of the second layer 
can continue in the same way as the first layer (Fig.[3l^f)- 
(g)). This contrasts with our previous model for vic- 
inal surface growth where a second layer of graphene 
can only grow at a step covered by a continous graphene 
strip (see Fig. [TJh)). When the propagation on the facet 
reaches the top junction between the nano-facet and the 
SiC(OOOl) basal plane, the graphene growth is allowed 
to continue on the (0001) basal plane but with a slower 
propagation rate Tp^^p = iyoex.p{-E^^^^/kT) (Fig. [S^d)- 
(e)). This is due to the experimental observation that 
several layers of graphene often grow on a nano-facet be- 
fore the graphene growth propagates on the (0001) basal 
plane^. We focus on the early growth stage where the 
graphene layers of each macrostep grow independently, 
and no coalescences or climb-overs are allowed. 

Fig.|3]shows that a particular choice of KMC model pa- 



3 



xlO-3 



AE/kT=0 



(a) 



' prop 




/I. /I 



AE/kT=5.8 



(b) 



xlO-3 


(c) 

AE/kT=7.7 








xlO-3 


(d) - 

AE/kT=11.6 









200 400 600 800 

Strip width, s 



1000 



FIG. 2. Graphene strip width distribution p{s) for different 
and total coverages with (j) — 0.9°. Different color lines 
correspond to G = 0.1 (red), 0.3 (magenta), 0.5 (green) and 
1.0 (blue), respectively. The terrace width is VK = 200. 



rameters produces a simulated morphology very similar 
to a transmission electron microscope image of graphene 
growth on a non-planar SiC(OOOl) surface. The TEM 
image was taken of a sample which was prepared at 1325 
°C for 90 min, with step heights of 5 — 15 nm. Both the 
TEM and the simulation image show a sharp ending for 
all graphene layers at the bottom of the nano- facet, where 
the graphene growth starts. Fig. |4] also shows that, from 
the left to the right side, graphene layers grow continu- 
ously over the junction between the high index nano-facet 
and SiC(OOOl). The top graphene layer is the longest. 
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FIG. 3. Graphene growth kinetics processes on a nano-facet. 
The blocks specified by letter A and B are used to calculate 
the fracture angle 0. 
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FIG. 4. (a) KMC simulation snapshot. Total surface coverage 
e = 0.8, T = 1800 K, Enuc/kT = 7.7, E'^,^^/kT is 3.9 and 
7.1 for the top graphene layer and the interface graphene lay- 
ers, respectively, (b) Transmission electron microscope image, 
adapted from Robinson et al.^. 



while the other layers have a similar length. In what fol- 
lows, we call all the graphene layers other than the top 
one "interface graphene layers" . 

Fig. [5] illustrates the formation process for the sur- 
face morphology in Fig. |4j The original nano-facet is 
shown in Fig. [5fa). The first graphene layer nucleates at 
the bottom junction of the nano-facet, and propagates 
quickly up the slope until it reaches the top basal plane. 
This is shown in Fig. [5jb). Fig. [5jc) shows that the first 
graphene layer continues to grow onto the terrace and a 
small triple bilayer nano-facet is fractured from the orig- 
inal nano-facet (This process was previously shown in 
Fig. [3l^d)-(e)). However, the terrace propagation rate is 
slower than the nano-facet propagation rate. Therefore, 
Fig. [5fc) also shows that the first graphene layer can- 
not grow very long on the terrace before the second layer 
is nucleated under the first layer and quickly covers the 
lower nano-facet. 

The terrace propagation rate for the second graphene 
layer is slower than the first layer due to the increasing 
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FIG. 5. Formation processes for the surface morphology in 
Fig.H 



difficulty for Si atoms to escape from SiC steps covered 
by graphene. This leads to a possible scenario, shown in 
Fig. Efd), that before the second graphene layer grows 
on the terrace, a third graphene layer is nucleated and 
catches up to the second layer growth at the top of the 
lower nano- facet. Later, when the second graphene layer 
does grow onto the terrace like the first layer, another 
triple bilayer step is decomposed from the lower nano- 
facet. This process is the same for the third layer. How- 
ever, assuming that the second and third graphene layers 
have a similar terrace propagation rate. Fig. [5fe) shows 
that a nano-facet of one unit cell high is fractured from 
the original nano-facet by growing the second and third 
graphene layers together on the terrace. 

Two small nano-facets are now fractured from the ori- 
gianl nano-facet, with heights corresponding to one and 
two triple bilayer steps, respectively. This is the same 
type of surface morphology seen in Fig. IH^b). We note 
that, in our model, graphene grows by decomposing SiC 
triple bilayer steps, so the graphene layer always has one 
side attached to a nano-facet. To emphasize this, in 
Fig. O we draw the graphene layers as slightly curved 
at the top of a nano-facet. 



IV. FRACTURE ANGLE DISTRIBUTION 

To quantify the surface morphology, we study the frac- 
ture angle made by the graphene layers on the nano- 
facet with respect to the basal plane. In the simulations, 
is quantified as follows. If the nano-facet height is /i, we 
count the number of units from the edge of the top terrace 
(marked as letter A in Fig.[3fg)) to the edge of the bottom 
terrace (letter B in Fig.[3l^g)), and define the length as L 
(e.g. L = 5 in Fig. Wig))- Then 6 = arctan((/i + 1)/L). 
The initial in the simulations is assumed to be Oq = 30°, 
which is close to the nano-facet angle observed experi- 
mentally. As soon as the top graphene layer growth prop- 



agates onto the terrace, starts to decrease as a func- 
tion of the growth time. Hence, we can use 6 to charac- 
terize the nano-facet fracturing process due to graphene 
growth. can be measured directly by experiments using 
the definition above. Compared to experimental results, 
the KMC simulation can be used to provide some quanti- 
tative information about the nucleation and propagation 
barriers. 

The initial surface for our statistical study of growth- 
induced fracturing is a periodic sequence of basal planes 
and nano-facets. We fix the height of the nano-facets 
to be 22.5 nm, which corresponds to = 30 SiC triple 
bilayers. We use a long terrace width and a short total 
growth time so that no graphene strip coalescence occurs. 
For statistical purposes, the surface consists of 6 x 10^ 
alternations of SiC(OOOl) and nano-facets, with a vicinal 
angle ^ = 5.7°. The growth temperature is fixed at 1800 
K. We also fix the nucleation barrier E^^c/kT = 7.7 and 
treat the terrace propagation barrier £^prop as the only 
variable. 

We now focus on the distribution of fracture angles at 
a given coverage. For simplicity only, we assume the en- 
ergy barriers for all graphene layers are the same. Fig. [6] 
shows the distribution p{0) of the normalized fracture an- 
gle tan (6)/ tan (^o) for different choices of terrace propa- 
gation barriers £^prop with a fixed total coverage 6 = 0.5. 
If we define AEf = E^mc — ^prop' clear from Fig. [6] 
that the distribution of fracture angles shows a cross-over 
behavior as A£^j increases. For A£^j = 0, it is difficult 
for graphene to propagate on the terrace. After all the 
nano-facets are fully covered by graphene, the growth on 
the terrace starts almost simultaneously. This leads to a 
Gaussian distribution of fracture angles in Fig. [6l This 
growth regime is very similar to graphene growth on vici- 
nal surfaces for /S.E = 0, where all the steps are nucleated 
at the same time and a Poisson strip width distribution 
is found. 

More interestingly, as AEf increases, Figl6] shows that 
the distribution of fracture angles becomes more and 
more one-sided. Eventually, when AEf = 5.8, the dis- 
tribution follows an intriguing power-law form, in which 
the probability is zero for values smaller than the min- 
imum 6 given in the figure. We do not fully under- 
stand the power-law probability distribution. Neverthe- 
less, the change in fracture angle distribution can still 
be related to the strip width distribution for vicinal sur- 
faces in Fig. O For growth on vicinal surfaces, when 
the nucleation barrier is high, graphene is nucleated at 
some steps much earlier than at the others, leading to 
a shift of peak position to the right side in the strip 
width distribution (Fig. Efb)). In other words, longer 
strips start to dominate the distribution. Similarly, for 
growth on nano-facets, when AEf > 0, the nano-facets 
are fractured one after another in a wide time range. In 
the regime where AEf is big enough, the distribution 
is dominated by smaller values of made by longer top 
graphene layers. This leads to a one-sided distribution of 
0. 
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FIG. 6. Fracture angle distribution for different terrace prop- 
agation barriers. The red dashed hne is a power- law fit, and 
the black dashed line is a Gaussian fit. 
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FIG. 7. Om as a function of KMC time. Dashed curves are 
the fit results according to Eq.[T] Red, green, blue and black 
dashed curves give E'-^^^^/kT — 2.1,4.0,6.0 and 8.0, respec- 
tively. 



Given the above discussion, similar to the strip width 
distribution, the cross-over behavior of the fracture angle 
distribution can be considered as a competition between 
the nucleation process at the nano-facet and the propa- 
gation process on the terrace. When compared to exper- 
imental observations, the fracture angle distribution can 
be used to determine A£^/. 

We define Om as the most probable in the fracture 
angle distribution. In the following, we give a simple 
analytic treatment to show how 9m evolves as a function 
of the growth time. If the experimental growth time is 



Ie-, the normalized 6m can be approximated by 
tan {Om) h^l 



tan ( 



(1) 



In this expression, h is the original nano-facet height, and 
^prop is effective terrace propagation barrier. To de- 
rive this equation, we do not consider the initial growth 
time spent on the nano-facet before the terrace growth 
occurs, which gives a small reduction to Ie in Eq. [H 
^prop accounts for both the top graphene layer terrace 
propagation barrier E'^^^^ and the nucleation barrier for 
interface graphene layers. Increasing the number of in- 
terface graphene layers will reduce the relative horizontal 
distance between the top and the bottom terraces at the 
fractured nano-facet (i.e. the horizontal distance between 
block A and B in FiglSfg)). Both the initial growth time 
spent on the nano-facet and the nucleation of interface 
graphene layers have the effect of increasing 0. There- 
fore, it is expected that ^p^op > ^prop- We treat ^p^op 
as a fitting parameter for our analytic model to see how 
it is related to ^prop- 

Fig. [71 shows the simulation result for the normal- 
ized Om as a function of the KMC time tk for differ- 
ent choices of terrace propagation barriers. Here we 
define A = E'^^^^/kT. The KMC time is defined as 
Tk — tE^prop' Dashed lines are fitted curves according to 
Eq. [TJ Despite the simphcity of Eq. [TJ the model curves 
agree quite well with the KMC simulations for smaller A 
and later growth times. This is likely caused by the ad- 
justment of initial growth time spent on the nano-facet. 
Smaller A indicates a faster graphene propagation on the 
terrace, which leads to a relatively shorter growth time 
spent on the nano-facet before the terrace growth occurs. 
Similarly, at later growth times when the initial growth 
time on the nano-facet becomes insignificant, the fitted 
curve also agrees better. 

As we increase the terrace propagation barrier ^p^op? 
the fitting parameter ^p^op is also found to increase. 
In fact, we find an excellent linear relationship between 
^prop and E'^^^^: E'^^^^ = E^^^^ + 0.03 eV. The rela- 
tively small correction to E'^^^^ suggests that the major 
contribution of £^prop still comes from £^prop foi" the top 
graphene layer. The contribution from the nucleation 
barrier for interface graphene layers acts as an additional 
energy barrier A£^prop = 0.03 eV. Therefore, with Eq. [H 
the time evolution of the fracture angle Om measured in 
an experiment can be used to obtain a good estimate of 
^prop- Combined with the result for AEf from the anal- 
ysis of fracture angle distribution, both £^nuc and E'^^^^ 
can be extracted experimentally. 

Most experiments of graphene growth on SiC nano- 
facet s are done at a temperature between 1600 K and 
1800 K. We tested A^^prop in this range, but no tem- 
perature dependence was found. A£^prop also does not 
seem to be dependent on the original nano-facet height 
h. As we increase the nano-facet height from 22.5 nm 
to 100 nm, the change of AEp^^p is less than 10%, well 
within the accuracy of £^prop as obtained from Fig. [3 
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This suggests that AEp^op is only a function of the en- 
ergy barriers in the model, and very weakly dependent on 
the growth temperature or the initial nano-facet height 
in the experimental parameter range. 



V. CONCLUSION 

We have presented a model for graphene growth on 
non-planar nano-facetted 6H-SiC substrates based on a 
previous model of growth on vicinal substrates. The sim- 
ulation produces a surface morphology very similar to 
observations. A description of the formation process for 
this type of surface morphology is also provided. For 
graphene growth on nano-facetted SiC substrates, a frac- 



ture angle can be used to characterize the growth-induced 
fracture of a nano-facet. The distribution of fracture an- 
gles is found to be related to the graphene strip width dis- 
tribution for growth on vicinal substrates. As the terrace 
propagation barrier decreases, the fracture angle distri- 
bution deviates from a Gaussian and eventually becomes 
a power-law. Our analytic result for the most probable 
fracture angle agrees well with the KMC simulations. 
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